# Comparing and synthesizing actor-elicited structures for complex food systems models
# Supporting code
# Feb 2021

library(dagitty)

a <- dagitty('dag {
  bb="-12.14,-9.869,13.385,9.559"
  Agri_inputs [pos="-2.398,-5.321"]
  Branding [pos="-1.440,1.730"]
  Brokers [pos="11.682,-8.890"]
  Cold_chain [pos="-3.037,4.781"]
  Cold_storage [pos="-6.231,6.566"]
  Consultants [pos="-3.703,-7.221"]
  Consumer [outcome,pos="10.377,-5.955"]
  Crop_mix [pos="-10.011,-6.099"]
  Demand [pos="9.818,5.961"]
  Depot [pos="5.294,-2.530"]
  Distribution [exposure,pos="-0.109,3.659"]
  Distribution_partnerships [pos="0.583,6.364"]
  Equipment [pos="-8.893,5.472"]
  Export [exposure,pos="10.191,-0.342"]
  Farm_labor [exposure,pos="-7.589,2.421"]
  Farms [pos="-4.741,-2.270"]
  Genetics [pos="-9.878,-8.257"]
  Geography [pos="-10.729,-1.637"]
  Governance [exposure,pos="0.210,-5.552"]
  Grading [exposure,pos="-1.493,-0.025"]
  Imports [pos="12.107,1.011"]
  Information [exposure,pos="-3.197,2.133"]
  Local_market [pos="11.469,-3.307"]
  Machinery [exposure,pos="-5.566,3.803"]
  Market_prices [pos="10.324,3.256"]
  Municipal_government [pos="6.012,-8.027"]
  National_government [pos="-0.003,-7.969"]
  Nursery [exposure,pos="-7.881,-8.516"]
  Packaging [exposure,pos="2.419,-1.954"]
  Planning [pos="-7.456,-0.486"]
  Politics [pos="8.594,-9.034"]
  Pre_negotiated_prices [pos="5.879,6.220"]
  Processing [pos="-2.026,-2.443"]
  Province_government [pos="3.085,-9.236"]
  Regulations [exposure,pos="2.845,-5.465"]
  Retail [pos="7.316,-4.918"]
  Shipping [pos="2.499,3.486"]
  Shipping_logistics [pos="5.267,3.975"]
  Soil [pos="-9.745,-3.335"]
  Storage [exposure,pos="2.872,1.011"]
  Supermarkets [exposure,pos="7.396,1.558"]
  Timing [pos="-10.012,0.887"]
  Training [exposure,pos="-10.516,3.457"]
  Water [exposure,pos="-5.646,-8.200"]
  Weather [pos="-2.558,-9.063"]
  Agri_inputs -> Farms
  Branding -> Supermarkets
  Brokers -> Consumer
  Brokers -> Local_market
  Brokers -> Retail
  Cold_chain -> Distribution
  Cold_chain <-> Cold_storage
  Consultants -> Farms
  Crop_mix -> Farms
  Demand -> Market_prices
  Depot -> Export
  Depot -> Retail
  Depot -> Shipping
  Depot -> Shipping_logistics
  Depot -> Supermarkets
  Distribution -> Depot
  Distribution -> Shipping
  Distribution_partnerships -> Distribution
  Distribution_partnerships -> Pre_negotiated_prices
  Equipment -> Cold_storage
  Equipment -> Machinery
  Export -> Local_market
  Farm_labor -> Farms
  Farm_labor -> Packaging
  Farms -> Grading
  Farms -> Processing
  Genetics -> Nursery
  Geography -> Crop_mix
  Geography -> Farms
  Geography -> Planning
  Geography -> Soil
  Governance -> Farm_labor
  Governance -> Regulations
  Grading -> Branding
  Grading -> Export
  Grading -> Local_market
  Grading -> Packaging
  Grading <-> Processing
  Imports -> Shipping
  Imports -> Supermarkets
  Information -> Planning
  Information -> Pre_negotiated_prices
  Local_market -> Consumer
  Local_market -> Supermarkets
  Machinery -> Farms
  Machinery -> Processing
  Market_prices -> Local_market
  Municipal_government -> Regulations
  National_government -> Regulations
  Nursery -> Farms
  Packaging -> Distribution
  Packaging -> Storage
  Planning -> Crop_mix
  Planning -> Farms
  Politics -> Export
  Politics -> Imports
  Politics -> Municipal_government
  Politics -> National_government
  Politics -> Province_government
  Pre_negotiated_prices -> Export
  Processing -> Depot
  Processing -> Packaging
  Province_government -> Regulations
  Regulations -> Distribution
  Regulations -> Farms
  Regulations -> Grading
  Regulations -> Packaging
  Regulations -> Retail
  Regulations -> Storage
  Retail -> Consumer
  Retail -> Supermarkets
  Shipping -> Export
  Shipping_logistics -> Shipping
  Soil -> Crop_mix
  Soil -> Planning
  Storage -> Retail
  Storage -> Shipping
  Supermarkets -> Consumer
  Timing -> Planning
  Training -> Farm_labor
  Water -> Crop_mix
  Water -> Farms
  Water -> Nursery
  Water -> Planning
  Water -> Processing
  Weather -> Water
}')

plot(a)

# Name all of the parents of Z.
parents( a, "Consumer" )

# Name all of the ancestors of Z.
ancestors( a, "Consumer" )

# Name all of the children of W.
children( a, "Farms" )

# Name all of the descendants of W.
descendants( a, "Farms" )

# Draw all (simple) paths between X and T.

paths( a, "Farms", "Consumer" )$paths

# Draw all the directed paths between X and T.

paths( a, "Farms", "Consumer", directed=TRUE )$paths

# For every variable V in the graph, find a minimal set of nodes that renders V independent of all other variables in the graph.

for( v in names( a ) ){
  cat( "Variable", v, ":", markovBlanket( a, v ), "\n" )
}

# Suppose we wish to estimate the value of "Consumer" from measurements taken on all other variables in the model. 
# Find the smallest set of variables that would yield as good an estimate as before.

predictors <- setdiff( names( a ), "Consumer" )
intersect( predictors, dconnected( a, "Consumer", list(), predictors ) )

# In general, this set would be the Markov blanket of Y. 
# In this specific case, the Markov blanket is just the set of the parents of Y.

# network analysis

library(igraph)

# create table of connections for each stakeholder
nodes <- read.csv("F4C stakeholder maps - Nodes.csv")
connections <- read.csv("F4C stakeholder maps - Connections.csv")
connections$ID1 <- paste0(connections$Stakeholder,connections$From) # create unique ID per participant
connections$ID2 <- paste0(connections$Stakeholder,connections$To) # create unique ID per participant

# create table of nodes for each stakeholder            
Community_garden_nodes <- nodes[nodes$Stakeholder=="Community_garden",
                                c("Node","Sticker_percent")]
Farmers_nodes <- nodes[nodes$Stakeholder=="Farmers",
                       c("Node","Sticker_percent")]
Informal_settlement_nodes <- nodes[nodes$Stakeholder=="Informal_settlement",
                                   c("Node","Sticker_percent")]
Local_gov_nodes <- nodes[nodes$Stakeholder=="Local_gov",
                         c("Node","Sticker_percent")]
Supermarket_nodes <- nodes[nodes$Stakeholder=="Supermarket",
                           c("Node","Sticker_percent")]
Ward_committee_nodes <- nodes[nodes$Stakeholder=="Ward_committee",
                              c("Node","Sticker_percent")]

# create connection datasets for each stakeholder
Community_garden_connections <- connections[connections$Stakeholder=="Community_garden",
                                            c("From","To")]
Farmers_connections <- connections[connections$Stakeholder=="Farmers",
                                   c("From","To")]
Informal_settlement_connections <- connections[connections$Stakeholder=="Informal_settlement",
                                               c("From","To")]
Local_gov_connections <- connections[connections$Stakeholder=="Local_gov",
                                     c("From","To")]
Supermarket_connections <- connections[connections$Stakeholder=="Supermarket",
                                       c("From","To")]
Ward_committee_connections <- connections[connections$Stakeholder=="Ward_committee",
                                          c("From","To")]

# create network and calculate centrality for each node
Community_garden_net <- graph_from_data_frame(d = Community_garden_connections,
                                  directed = TRUE,
                                  vertices = Community_garden_nodes)
Farmers_net <- graph_from_data_frame(d = Farmers_connections,
                                              directed = TRUE,
                                     vertices = Farmers_nodes)
Informal_settlement_net <- graph_from_data_frame(d = Informal_settlement_connections,
                                              directed = TRUE,
                                              vertices = Informal_settlement_nodes)
Local_gov_net <- graph_from_data_frame(d = Local_gov_connections,
                                              directed = TRUE,
                                       vertices = Local_gov_nodes)
Supermarket_net <- graph_from_data_frame(d = Supermarket_connections,
                                              directed = TRUE,
                                         vertices = Supermarket_nodes)
Ward_committee_net <- graph_from_data_frame(d = Ward_committee_connections,
                                              directed = TRUE,
                                            vertices = Ward_committee_nodes)


# compare networks

tgt <- Informal_settlement_net # set target network

c1 <- eigen_centrality(tgt,
                       directed = FALSE,
                       scale = TRUE)

length(c1$vector) # N nodes
c1$vector # centrality for each node
sort(c1$vector) # centrality for each node in ascending order
mean(c1$vector) 
sd(c1$vector)

d1 <- degree(tgt, 
               mode = "all")

d1 # degree for each node
sort(d1)

edge_density(tgt, loops = TRUE)
transitivity(tgt, type = "global")
c1 <- cluster_optimal(tgt)
c1[1:length(c1)]

a1 <- Informal_settlement_nodes[order(Informal_settlement_nodes$Sticker_percent),]
a1
#


# visualise synthesised networks
#synth_nodes <- nodes[, c("Node","Sticker_percent")] # gives duplicates
synth_connections <- connections[,c("From","To")]

synth_net <- graph_from_data_frame(d = synth_connections,
                                              directed = TRUE)
#vertices = synth_nodes)

net <- simplify(synth_net)

#V(net)$color <- colrs
#deg <- degree(net, mode = "all")
com <- cluster_spinglass(net, spins=2)
V(net)$color <- com$membership+1
#V(synth_net)$size <- net$Sticker_percent
V(net)$label <- NA
#E(net)$width <- E(net)$weight/6

plot(net,
     layout=layout_nicely)
